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Abstract 

Molecular dynamics (MD) methods are opening new opportunities for simulating the 
fundamental processes of material behavior at the atomistic level. However, current analysis is 
limited to small domains and increasing the size of the MD domain quickly presents intractable 
computational demands. A preferred approach to surmount this computational limitation has 
been to combine continuum mechanics-based modeling procedures, such as the finite element 
method (FEM), with MD analyses thereby reducing the region of atomic scale refinement. Such 
multiscale modeling strategies can be divided into two broad classifications: concurrent 
multiscale methods that directly incorporate an atomistic domain within a continuum domain and 
sequential multiscale methods that extract an averaged response from the atomistic simulation 
for later use as a constitutive model in a continuum analysis. 

In the present work, new approaches to both concurrent and sequential multiscale modeling 
are developed. A new approach to concurrent MD-FEM coupling uses statistical averaging of 
the atomistic MD domain to provide displacement interface boundary conditions to the 
surrounding continuum FEM domain, which, in return, generates interface reaction forces 
applied as piecewise constant traction boundary conditions to the MD region. Additionally, a 
new approach to sequential multiscale modeling is developed and based on the definition of a 
cohesive zone volume element (CZVE). The CZVE is an atomistic analog to the cohesive zone 
models that are often used to model material degradation in continuum finite element analyses. 
Taken together, these two approaches are an attempt to recast atomistic processes within a 
format suitable for inclusion in grain-scale finite element analyses and are providing a necessary 
step toward true physics-based analysis of damage processes. 


I. Introduction 

Classical fracture mechanics is based on the comparison of computed fracture parameters to 
their empirically-determined critical values. Although this paradigm has been extremely 
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successful for modeling crack growth at structural scales, it does not describe the physical 
processes that govern fracture. Ultimately, the characterization of events and processes that 
occur at length scales on the order of 10~ 9 to 10" 3 m is needed to understand the physics of crack 
growth. 


Physics-based modeling of fracture begins at nanometer length scales in which atomistic 
simulation is used to predict the formation, propagation, and interaction of fundamental damage 
mechanisms. These mechanisms include dislocation formation and interaction, interstitial void 
formation, and atomic diffusion. These nanoscale damage mechanisms progress into microscale 
processes such as local plasticity and small crack formation. Ultimately, damage progression 
leads to macroscopic failure. 

However, modeling atomistic processes quickly becomes computationally intractable as the 
system size increases. With current computer technology, the computational demands of 
modeling suitable domain sizes (on the order of hundreds of atoms for quantum mechanics-based 
methods, and potentially billions of atoms for classical mechanics-based methods) and 
integrating the governing equations of state over sufficiently long time intervals quickly reaches 
an upper bound for practical analyses. In contrast, continuum mechanics methods such as the 
finite element method (FEM) provide an economical numerical representation of material 
behavior at length scales in which continuum assumptions are appropriate. Thus, multiscale 
modeling strategies are needed that use the most efficient material model at each length scale to 
develop a unified description of the hierarchy of processes that govern fracture. 

Multiscale modeling strategies relevant to length scales between the atomic-scale and the 
lower end of the grain-scale can be divided into two broad classifications: concurrent multiscale 
and sequential multiscale. Concurrent multiscale modeling bridges length scales by coupling 
different computational paradigms and is attractive as a highly efficient means of reducing the 
computational cost of the simulations in cases that require modeling of relatively large material 
domains to capture the complete deformation field, but where atomic and subatomic refinement 
is needed only in very localized regions (e.g., near a crack tip or dislocation core). Such 
computational issues arise in modeling crack nucleation and propagation, and in modeling 
dislocation formation and interaction. Concurrent multiscale models allow the size limitations of 
the atomistic simulation to be avoided by embedding an inner MD region (where complex 
dynamic processes and large deformation gradients exist) within an outer domain (where small 
deformation gradients exist) where a continuum representation of the material becomes 
appropriate. 

In contrast, sequential multiscale modeling relies on recasting previously developed 
characterizations of damage processes as constitutive models within higher-level analyses. 
Cohesive zone models have often been used as a means of sequential multiscale modeling by 
developing an interfacial traction-displacement relationship using bulk material properties to 
model possible fracture surfaces in a continuum. 1 " 5 The use of cohesive zone models in finite 
elements permits highly discontinuous displacement behavior to be exhibited in a single 
discretized model without remeshing while using an incremental piecewise linear solution 
algorithm. 
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In this paper, elements of both concurrent and sequential multiscale modeling will be used to 
illustrate the possibilities of characterizing fracture along grain boundaries in polycrystalline 
aluminum. An outline of a new concurrent multiscale method for coupling MD and FEM 
computational domains is presented in Section II. Section III summarizes the determination of 
cohesive zone relationships developed from molecular dynamics analysis using a new sequential 
multiscale modeling method. Also in Section III is a discussion of a simulation of an edge crack 
propagating through a continuum into an atomistic domain with the deformation along the crack 
plane in the continuum region modeled using CZMs. Section IV presents a simulation of 
intergranular fracture in a polycrystal where the CZMs are used as constitutive models along 
idealized grain boundaries. Finally, Section V presents a summary of the analysis and overall 
methodology. 


II. Concurrent Multiscale Modeling 

Over the past decade, various methods that directly couple material representations at 
different levels of refinement have been developed and offer significant computational 
advantages compared to full atomistic simulations for predicting deformation and fracture 
processes. 6 ' 16 A common feature of all of these approaches for coupling atomistic regions with 
continuum domains is the refinement of the FEM mesh to the atomic scale. While such direct 
coupling approaches are conceptually straightforward, the fundamental difficulty in their 
development lies in the inherent differences between the atomistic and continuum computational 
models. The physical state of the atomistic region is described through nonlocal interatomic 
forces between discrete atoms of given position and momentum, while the physical state of the 
continuum region is described through continuous stress-strain fields that reflect local statistical 
averages of atomic interactions at larger length and time scales. 

However, a robust connection between continuum and discrete quantities can be achieved 
through a statistical averaging over length and time scales where the discreteness of the atomic 
structure can be successfully approximated as a continuum. With the use of statistical averages 
to couple the two computational schemes, the developed approach constitutes a statistical 
coupling approach, and the developed MD-FEM coupling method is referred to as the embedded 
statistical coupling method (ESCM). 17 ’ 18 The resulting model that consists of an MD system 
embedded in a FEM domain is depicted in Figure 1 . 

The structure of the basic ESCM model consists of four regions that are depicted in Figure 2. 
The regions include the inner MD region, interface MD/FEM region, surface MD region and 
FEM domain. The inn er MD region is used to model the damage processes (e.g., plastic 
deformation events such as dislocation or void nucleation, crack propagation) that evolve during 
the simulation. The interface MD/FEM region provides the computational linkage between 
statistical averages of atomistic quantities from the MD and continuum forces and displacements 
from the FEM. The surface MD region does not interact with superposed FEM nodes but, rather, 
is used to compensate for atomic free edge effects at the outer boundary of the MD region. 
Finally, the FEM domain permits a large reduction in the computational cost of the simulations 
by replacing the atomistic representation with a continuum model in those regions of the system 
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where the deformation gradients are small and atomic-level resolution is not necessary. 
Complete details of the basic method are presented by Saether et al. 17, 18 

The basic algorithm of ESCM involves the statistical averaging over both time and volume of 
atomistic subdomains at the MD-FEM interface that provide nodal displacement boundary 
conditions to the continuum FEM domain. An analysis of the FEM domain is performed, which, 
in turn, generates interface reaction forces that are applied as constant traction boundary 
conditions 19 to the atoms within the localized MD subdomain. The MD system is simulated to 
accommodate the interface reaction forces and to provide an update of the displacement 
boundary conditions at the MD-FEM interface. This process is iterated until equilibrium between 
the atomistic MD forces and the FEM reaction forces is established at the interface. Typically, 
one finite element at the interface encompasses a region of several hundred to several thousand 
atoms. At this scale, the discreteness of the atomic structure is homogenized enough so that the 
FEM domain responds to the atomistic domain as an extension of the continuum. 

A recent addition to the ESCM involves the addition of cohesive zone model (CZM) elements 
along the crack plane to improve the compatibility of the deformations near the MD to FEM 
interface. Cohesive zone models approximate debonding processes along an interface 1 ' 5, 20 and 
are frequently used in conjunction with finite element approaches to study fracture in a wide 
variety of materials. Recently, cohesive zone volume elements (CZVEs) 21 have been introduced 
in the MD analysis to calculate atomistic-based parameters defining CZMs. The resulting CZMs 
exhibit the aggregate behavior of a volume of atoms in the MD simulation. 21 By placing the 
CZM elements along the plane of crack growth within the finite element domain of the ESCM 
model, the fracture process can be simulated in the continuum as the crack approaches the MD- 
FE interface. Figure 3 depicts the coupled MD-FEM model containing an edge crack. CZM 
elements are identified that represent the cohesive properties of the grain boundary (GB) in the 
continuum. Crystal I and crystal II identify the MD regions that have been rotated through a tilt 
angle to form the GB. The CZVEs will be discussed in more detail in the Section III. 

The simulation presented in this paper considers intergranular fracture wherein MD 
simulations of fracture along grain boundaries (GBs) are used to obtain the atomistic response of 
crack propagation under applied loads along a characteristic length of the GB. Intergranular, 
rather than transgranular, fracture is considered because the crack path is well-defined and is 
known a-priori. 

Figure 4 shows a <11 0> 299 symmetric tilt GB, a typical high energy GB formed between 
two grains in a fee metal. Because of the slip plane orientation near this GB, crack growth in the 
+x direction was observed to be brittle with few dislocations emitted, while crack growth in the - 
x direction was observed to be ductile with many dislocations and twins emitted. 21 Common 
neighbor analysis (CNA) 22 ’ 23 is used to identify atoms in different crystallographic states: fee 
(small dots), hep (triangles), non-crystalline atoms (large dots), and surface atoms (squares), 
indicating existing vacancies in the GB. The length scale is in units of the lattice constant of Al, 
a 0 = 0.405 nm, and the analysis is performed using the potential of Mishin and Farkas. 24 
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III. Sequential Multiscale Modeling 


Development of Traction-Displacement Relationships from MD 

Through the definition of CZVEs shown in Figure 5, the results from MD simulations have 
been recast to obtain averaged continuum traction-displacement relationships that represent the 
cohesive interactions along a characteristic length of material. Using a purely MD model of a 
center crack along the 299 grain boundary presented in Reference 21, the stress and opening 
profiles are extracted along the crack growing in the system as shown in Figure 6. The figure 
shows distinct traction, cf yy (x), and displacement, A, profiles for the two crack tips of the center 
crack. The asymmetric behavior results from the orientation of the crystal lattice of the grains 
forming the grain boundary. A curve fit of x m in Figure 6 is shown as a close approximation to 
the cfyy{x) values extracted from the MD for the brittle tip. In contrast, plastic deformation 
processes, including twinning and dislocation emission, were found to dominate deformation at 
the ductile crack tip growing in the opposite direction. 21 

The cfyy(x) and A(x) profiles in Figure 6 represent the traction and opening displacement at a 
single instant of time during crack propagation. If the entire simulation of the interface 
debonding is divided into N t equal intervals of time t q (<?E[0, N t ]), numerous similar profiles can 
be taken of many CZVEs placed along the GB. When plotted as a d yy vs. A curve, each 
0 cfyyixpjq ), Mx p , t q )) pair represents a point d vv (X(xp,t q ))- After sorting these data points in order 
of ascending A, so that A, < A ,•+/, and taking a moving average (or a consecutive mean), 


t (g)= 


(2M + 1), 


M / \ 


j—M 
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in which the results are averaged over M points backward and M points forward from A„ a 
construction of a statistically representative cohesive zone relationship t(A) can be made. 
Complete details of the modeling, fracture processes, and extraction of the CZM response are 
given in Yamakov, et al. 21 

Crack Propagation through a Continuum-Atomistic Interface 

A continuum-atomistic model for crack propagation in the brittle direction along the 299 
GB shown in Figure 3, and analyzed using the enhanced ESCM analysis, is presented. The 
dimensions of the coupled MD-FEM model used for this simulation are r/pp = 900 nm in the x- 
and y- directions, with a circular MD region in the x-y plane of diameter, Jmd = 45 nm, as shown 
in Figure 1. The system thickness in the z-direction is h = 2.9 nm with periodic boundary 
conditions to emulate a bulk atomic state in this direction. The simulation was performed at a 
temperature T= 300 K and maintained by a Nose-Hoover thermostat. 25 

An edge crack was inserted along the line of symmetry, c, in the FEM mesh which ended 
with four finite elements on the left side before reaching the MD-FEM interface (Figure 3). 
These four remaining FEs were connected through CZM elements placed along the crack 
direction (line c in Figure 3), which meet the atomistic GB at the MD-FEM interface. On the 
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right side of the MD region, the debonding interface is continued through another line of CZM 
elements, which extends to the opposite outer boundary of the FEM domain. 


The system was subjected to uniaxial far-field strain, %, applied as displacement boundary 
conditions at the outer boundary of the FEM domain. The strain was gradually increased from 
0.4% to 2% in steps of 0.4% every 20 MD-FEM iterations (20 ps MD simulation time), allowing 
the MD system to more smoothly follow the stress increase in the surrounding FEM domain. The 
crack grows as an edge crack in the continuum environment subjected to mode I loading 
conditions and penetrates into the atomistic environment (Figure 7). The continuity of the crack 
propagation process through the MD-FEM interface is achieved by using the CZM constitutive 
relation derived for interface debonding based on the pure molecular-dynamics simulations 21 as 
discussed previously. Figure 7 presents the system configuration and o yy stress state at the end of 
the second iteration when £ yy = 2.0%. 

The process of crack propagation is illustrated in Figure 8 through a set of snapshots, 
magnified to show the physical configuration of the MD region in the center of the MD-FEM 
coupled model. The GB interface in the MD region is visualized as a line of blue and red atoms, 
indicating atoms in a disordered or hep surrounding, respectively, reflecting the non-crystalline 
structure of the GB interface. In Figures 8(a) and 8(b), the crack propagates through the FEM 
mesh, sequentially opening the CZM elements placed on its path. In Figure 8(c), the crack tip 
enters the MD domain. The snapshot has captured the instant when the debonding of the surface 
and interface MD Regions has begun in unison with the opening of the last CZM element on the 
left side of the MD-FEM interface. The crack tip is about to enter the inner MD region and to 
continue its propagation in a fully atomistic way. In Figure 8(d), decohesion of the GB in the 
inner MD region has begun. Shortly afterward, in Figure 8(e), plastic processes start to develop 
in terms of nucleation of a pair of dislocations, seen as two short blue lines, marking the 
extended dislocation cores on both sides of the crack tip. Figure 8(f) presents the crack 
propagated well inside the MD domain. Active dislocation nucleation has already developed and 
accompanies the crack’s propagation. An important result of this process, shown in Figures 8(a- 
f), is the preservation of the structural continuity between the atomistic and the continuum 
domains at the MD-FE interface through the entire process of crack propagation from the FEM 
to the MD systems. 

The issues related to coupling the dynamic response of the MD domain to the static response 
of the FEM system are discussed in detail in Saether, et al. 17,18 While the dynamic coupling in the 
ESCM approach cannot be fully achieved without the use of a dynamic FEM simulation, the 
iterative scheme of coupling, in which the FEM state is continuously updated in accordance with 
the reported MD displacements, appears to be sufficient to reproduce the evolution of systems 
with relatively slow dynamics. 

Development of Traction-Displacement Relationships from ESCM 

In the original purely MD simulation 21 that was used to obtain CZM relationships, model size 
restrictions resulted in a finite geometry distortion. The model was limited to a domain of 
approximately 120 nm x 97 nm x 3 nm and required the application of extremely large applied 
stresses (approximately 4 GPa) to grow an initial crack that was only a few nanometers in length. 


Page 6 of 17 

Glaessgen, et al., 2008 Aging Aircraft Conference 



The hydrostatic stress was applied corresponding to fixed displacement boundary conditions 
before the crack was introduced. Because of the extremely large stresses, hydrostatic loading 
was required to suppress bulk plasticity. 

In contrast, the enhanced ESCM simulation allows uniaxial far-field loading to be applied 
and eliminates the finite geometry distortion but has only been implemented for the “brittle” 
crack tip. The procedure for the determination of the CZM parameters of the brittle crack tip 
using the coupled MD-FEM simulation is the same as that used in the purely MD simulation. 
The insert in the bottom half of Figure 9 shows a schematic diagram of the discretization of the 
system volume in the MD simulation into CZVEs for extracting the parameters for the CZM 
elements. As in Figure 6, Figure 9 shows the cf yy (x) and A(x) profiles that represent the traction 
and opening displacement at a single instant of time during crack propagation. The stress profile 
immediately ahead of the crack tip follows a 1/Vx type distribution. 

A comparison of the recomputed brittle Mode I CZM with the initial pure MD results from 
Reference 21 is presented in Figure 10. The ESCM-determined traction-displacement curve 
shows a very large but finite slope until the peak traction of approximately 6 GPa is reached. In 
contrast, the MD-determined curve begins at approximately 4.25 GPa, corresponding to the 
value of the hydrostatic pre-stress used in the simulation, and reaches a peak value of 
approximately 5 GPa. The remainders of both curves consist of a softening region followed by 
complete unloading (corresponding to the formation of a stress-free surface) at just over 2 nm. 


IV. Grain-Scale Idealization 

The CZM relationships extracted from the ESCM analysis contain many idealizations 
(e.g., absence of defects) along the crack plane and cannot be deemed as being quantitative; 
however, they do improve our understanding of the underlying mechanisms of damage at very 
small length scales. Additionally, since the selected Z99 GB is believed to be representative of a 
common class of GBs, additional insight into the mechanisms of intergranular fracture may be 
gained by implementing the CZMs within a microstructural analysis. 

Figure 1 la shows a demonstration problem consisting of a 2-D configuration of a tessellated 
polycrystalline microstructure having many grain boundary interfaces. The model dimensions 
are l x = l y = 47.5 nm encompassing 50 grains. 26 CZM parameters characterizing the brittle Mode 
I crack tip behavior of the MD-generated results were extracted from the relationships shown in 
Figure 10 and (for the purposes of this demonstration) are assumed to apply to all grain 
boundaries regardless of orientation. These parameters are shown in Table 1, and the 
constitutive behavior of the resulting idealized CZM is shown by the red curve in Figure 1 lb. 


Table 1 . CZM parameters approximating brittle Mode I behavior at the crack tip 


CZM 

Brittle Mode I 

6"„,v (GPa) 

5.00 

6° (nm) 

0.2 

S mx (nm) 

2.17 
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Using brittle Mode I CZM properties for all the GBs, Figure 12 shows the sequence of 
decohesion processes under normal applied displacements with the opening being exaggerated in 
the figure. An examination of the decohesion processes exhibited in the simulation presented in 
Figure 12 shows the initial opening of GBs under normal applied displacements. The initial 
opening of the GBs is shown in Figure 12a. At higher applied strains, the energetics of GB 
opening dictate the location at which a dominant crack begins to form (Figure 12b). At a critical 
applied strain, the local coalescence of GB opening forms a dominant microcrack as shown in 
Figure 12c. After the dominant microcrack has fractured the model, the elastic stiffness in other 
CZM elements acts to return them to an unloaded state. 


V. Summary 

New approaches to both concurrent and sequential multiscale modeling have been developed. 
A new approach to concurrent MD-FEM coupling uses statistical averaging of the atomistic MD 
domain to provide displacement interface boundary conditions to the surrounding continuum 
FEM domain, which, in return, generates interface reaction forces applied as piecewise constant 
traction boundary conditions to the MD domain. Additionally, a new approach to sequential 
multiscale modeling was developed and based on the definition of a cohesive zone volume 
element (CZVE). The CZVE is an atomistic analog to the cohesive zone models (CZM) that are 
often used to model material degradation in continuum finite element analyses. Taken together, 
these two approaches are an attempt to recast atomistic processes within a format suitable for 
inclusion in grain-scale finite element analyses and are providing a necessary step toward true 
physics-based analysis of damage processes. 

An application of the ESCM, enhanced by the addition of CZM elements along the crack 
plane, was presented for analysis of a Z99 GB in aluminum. The model incorporated an edge 
crack with the tip propagating into the inner MD region of ESCM. The continuity of the crack 
propagation process through the MD-FEM interface was achieved by the derivation of the CZM 
constitutive relation for interface debonding using CZVEs and moving averages based on an 
earlier pure molecular-dynamics simulation. CZM parameters for Mode I fracture were 
determined and compared with those determined from a pure MD simulation. For a “brittle” 
crack tip, the peak stress was increased slightly in the ESCM simulation but the overall energy 
and functional form of the CZM was closely reproduced. Finally, the MD-derived CZMs were 
used to model intergranular fracture of a highly idealized microstructure. 


VI. Acknowledgments 

V. Yamakov is sponsored through cooperative agreement NCC- 1-02043 with the National 
Institute of Aerospace. 


Page 8 of 17 

Glaessgen, et al., 2008 Aging Aircraft Conference 



VII. References 

'Klein, P. and Gao, H., “Crack Nucleation and Growth as Strain Localization in a Virtual- 
Bond Continuum f Engineering Fracture Mechanics, Vol. 61, 1998, pp. 21-48. 

2 Camacho, G.T. and Ortiz, M., “Computational Modeling of Impact Damage in Brittle 
Materials,” Int. J. Solids Struct., Vol. 33, 1996, pp. 2899-2938. 

3 Tvergaard, V. and Hutchinson, J.W., “The Relation Between Crack Growth Resistance and 
Fracture Process Parameters in Elastic-Plastic Solids,” J. Mech. and Phys. Solids, Vol. 40, 1992, 
pp. 1377-1397. 

4 Iesulauro, E., Ingraffea, A. R., Arwade, S. R., and Wawrzynek, P. A., “Simulation of Grain 
Boundary Decohesion and Crack Propagation in Aluminum Microstructure Models,” Fatigue 
and Fracture Mechanics: 33 rd Volume, ASTM STP 1417, 2001, pp. 715-728. 

5 Tvergaard, V. and Hutchinson, J.W., “The Influence of Plasticity on Mixed-Mode Interface 
Toughness,” J. Mech. Phys. Solids, Vol. 41, 1993, pp. 1119-1135. 

6 Abraham, F.F., Broughton, J., Bernstein, N., and Kaxiras, E., “Spanning the Continuum to 
Quantum Length Scales in a Dynamics Simulation of Brittle Fracture,” Europhys. Lett., Vol. 44, 
1998, pp. 783-787. 

7 Broughton, J., Abraham, F.F., Bernstein, N., and Kaxiras, E., “Concurrent Coupling of 
Length Scales: Methodology and Application,” Phys. Rev. B, Vol. 60, 1999, pp. 2391-2403. 

8 Curtin, W.A. and Miller, R.E., “Atomistic/Continuum Coupling in Computational Materials 
Science,” Modelling Simul. Mater. Sci. Eng., Vol. 11, 2003, pp. R33 - R68. 

9 Park, H.S. and Liu, W.K., “An Introduction and Tutorial on Multiple-scale Analysis in 
Solids,” Comput. Methods inAppl. Mech. Engrg., Vol. 193, 2004, pp. 1733-1772. 

10 Li, X. and E, W., “Multiscale Modeling of the Dynamics of Solids at Finite Temperature,” 
J. Mech. Phys. Solids, Vol. 53, 2005, pp. 1650-1685. 

"Miller, R.R. and Tadmor, E.B., “The Quasicontinuum Method: Overview, Applications and 
Current Directions,” J. Comp. Aided Mat. Desgn, Vol. 9, 2002, pp. 203-239. 

12 Xiao, S.P. and Belytschko, T., “A Bridging Domain Method for Coupling Continua with 
Molecular Dynamics,” Comput. Methods Appl. Mech. Engrg., Vol. 193, 2004, pp. 1645-1669. 

13 Shilkrot, L.E., Miller, R.E., and Curtin, W.A., “Coupled Atomistic and Discrete Dislocation 
Plasticity,” Phys. Rev. Lett., Vol. 89, 2002, p. 025501-1. 

14 Gumbsch, P. and Beltz, G.E., “On the Continuum Versus Atomistic Description of 
Dislocation Nucleation and Cleavage in Nickel,” Mod. Sim. Mat. Sci. and Engrg., Vol. 3, 1995, 
pp. 597-613. 

15 Rudd, R.E. and Broughton, J.Q., “Coarse-Grained Molecular Dynamics and the Atomic 
Limit of Finite Elements,” Phys. Rev. B, Vol. 58, 2005, pp. R5893-R5896. 

16 Rudd, R.E. and Broughton, J.Q., “Coarse-Grained Molecular Dynamics: Nonlinear Finite 
Elements and Finite Temperature,” Phys. Rev. B, Vol. 72, 2005, pp. 144104-1-32. 

17 Saether, E., Yamakov, V., and Glaessgen, E.H, “A Statistical Approach for the Concurrent 
Coupling of Molecular Dynamics and Finite Element Methods,” Proceedings of the 48th 
AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics, and Materials Conference, 
AIAA-2007-2 1 69-CP, AIAA, Honolulu, HI, April 23-26, 2007. 

18 Saether, E., Yamakov, V., and Glaessgen, E.H., “An Embedded Statistical Method for 
Coupling Molecular Dynamics and Finite Element Analyses,” 2008, (submitted to Int. J. Num. 
Meth. Engr.). 

19 Cleri, F., “Representation of Mechanical Loads in Molecular Dynamics Simulations,” Phys. 
Rev. B, Vol. 65, 2001, pp. 014107-1-6. 


Page 9 of 17 

Glaessgen, et al., 2008 Aging Aircraft Conference 



20 Costanzo, F. Allen, D.H., A Continuum Thermodynamic Analysis of Cohesive Zone 
Models,” Int. J. Sci. Eng. , Vol. 33, pp. 2 1 97-22 19. 

21 Yamakov, V., Saether, E., Phillips, D. R., and Glaessgen, E. H., “Molecular-Dynamics 
Simulation-Based Cohesive Zone Representation of Intergranular Fracture Processes in 
Aluminum,” J. Mech. Phys. Solids, Vol. 54, 2006, pp. 1899-1928. 

22 Honeycutt, J. D. and Andersen, H. C., “Molecular Dynamics Study of Melting and Freezing 
of Small Lennard-Jones Clusters,”/. Phys. Chem., Vol. 91, 1987, pp. 4950-4963. 

23 Clarke, A.S. and Jonsson, H., “Structural Changes Accompanying Densification of Random 
Hard-Sphere Packings,” Phys. Rev. E, Vol. 47, 1993, pp. 3975-3984. 

24 Mishin, Y. and Farkas, D., “Interatomic Potentials for Monoatomic Metals from 
Experimental Data and Ab Initio Calculations,” Phys. Rev. B, Vol. 59, 1999, pp. 3393 - 3407. 

25 Nose S., “A Unified Formulation of the Constant Temperature Molecular Dynamics 
Method,” J. Chem. Phys. , Vol. 8 1 , 1 984, p. 5 1 1 . 

26 Glaessgen, E.H., Saether, E., Phillips, D.R., and Yamakov, V., “Multiscale Modeling of 
Grain-Boundary Fracture: Cohesive Zone Models Parameterized from Atomistic Simulations,” 
Proceedings of the 47th AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics, and 
Materials Conference, AIAA-2006-1674-CP, AIAA, Newport, RI, May 1-4, 2006. 


Page 10 of 17 

Glaessgen, et al., 2008 Aging Aircraft Conference 



VIII. Figures 



Figure 1 . An MD region embedded within a FEM domain 



Figure 2. Structure of embedded statistical coupling method 
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Figure 3. Enhanced ESCM for edge-crack propagation along a grain boundary in 
aluminum 


[7 7 10] 



Figure 4. Crystallography and structure of the simulated 
GB interface in aluminum 
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Figure 5. Schematic diagram of the volume slicing and CZVE regions in the 
MD simulation of a 299 grain boundary 



Figure 6. Stress and opening profiles extracted along an embedded crack growing in a 
MD simulation. A snapshot of the corresponding crack configuration is shown at the 
bottom of the figure. 
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Figure 7. ESCM simulation of an edge crack propagating along a <1 10> 299 symmetric tilt grain 
boundary in aluminum for 2.0% far-field uniaxial strain 
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MD atoms color scheme: 

o fee; o MD interface fee; • hep; o disordered; o surface 


Figure 8. Snapshots of the ESCM model configurations simulating an edge crack propagating 
along a GB in a bicrystal of aluminum with CZM elements along the GB interface 
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Figure 9. Stress and opening profiles extracted from cohesive-zone-volume-elements 
(CZVE) placed along the path of an edge crack growing along a grain-boundary inside 
the molecular dynamics region 



Figure 10. Comparison between the brittle CZM derived using 
MD and ESCM models 
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(a) Idealized microstructure 


(b) Idealized CZM 


Figure 1 1 . Idealizations of microstructure and cohesive zone model 




(a) Initial opening of GBs (b) Local crack formation (c) Dominant microcrack 

£yy= 1 6% £yy =17% £yy = 22% 

Figure 12. Process of dominant microcrack formation in a polycrystal model 
simulated using idealized cohesive zone models 


Page 17 of 17 

Glaessgen, et al., 2008 Aging Aircraft Conference 





